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Abstract 

Various bacterial strains exhibit colonial branching patterns during growth on 
poor substrates. These patterns reflect bacterial cooperative self-organization 
and cybernetic processes of communication, regulation and control employed 
during colonial development. One method of modeling is the continuous, 
or coupled reaction-diffusion approach, in which continuous time evolution 
equations describe the bacterial density and the concentration of the relevant 
chemical fields. In the context of branching growth, this idea has been pursued 
by a number of groups. We present an additional model which includes a 
lubrication fluid excreted by the bacteria. We also add flelds of chemotactic 
agents to the other models. We then present a critique of this whole enterprise 
with focus on the models' potential for revealing new biological features. 



I. INTRODUCTION 



During the course of evolution, bacteria have developed sophisticated cooperative behav- 
ior and intricate communication capabilities [||-|^. These include: direct cell-cell physical 
interactions via extra-membrane polymers collective production of extracellular "wet- 
ting" fluid for movement on hard surfaces , long range chemical signaling, such as quorum 



1 



sensing [p!0|-[T2[] and chemotactic signaling^ [1T3|-[T5[], collective activation and deactivation of 
genes 



r^-p!8| and even exchange of genetic material [|T^-pT|. Utilizing these capabilities, 



bacterial colonies develop complex spatio-temporal patterns in response to adverse growth 
conditions. 

It is now understood that the study of cooperative self-organization of bacterial colonies 
is an exciting new multidisciplinary field of research, necessitating the merger of biological 
information with the physics of non-equilibrium processes and the mathematics of non- 
linear dynamics. At this stage, several experimental systems have been identified, and 
preliminary modeling efforts are making significant progress in providing a framework for 



the understanding of experimental observations ||8|JT8, 22-32 



Fujikawa and Matsushita p3| , |33| , p4| reported for the first time f\ that bacterial colonies 
could grow elaborate branching patterns of the type known from the study of fractal forma- 
tion in the process of diffusion-limited-aggregation (DLA) ||37|^9| . This work was done with 
Bacillus subtilis, but was subsequently extended to other bacterial species such as Serra- 
tia marcescens and Salmonella anatum . It was shown explicitly that nutrient diffusion 
was the relevant dynamics responsible for the growth instability. Later, we will see how 
models which couple nutrient diffusion to bacterial density can naturally account for these 
structures. 



Motivated by these observations, Ben- Jacob et al. ||25| , ^ , |27[| conducted new experiments 
to see how adaptive bacterial colonies could be in the presence of external "pressure" , here 
in the form of a limited nutrient supply and hard surface. The endeavor started with B. 
subtilis 168, which is non-motile on a solid agar surface, from which a new species of bacteria 



Chemotaxis is a bias of movement according to the gradient of a chemical agent. Chemotactic 
signaling is a chemotactic response to an agent emitted by the bacteria. 
^ We refer to the first time that branching growth was studied as such. Observations of branching 
colonies occurred long ago 
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has been isolated [^,0 • The new species was designated Paenihacillus dendritiformis var. 
dendron |^^. This species is motile on the hard surface and its colonies exhibit branching 
patterns (Fig. 1). The new mode of tip-splitting growth was found to be inheritable and 
transferable by a single cell, hence it is referred to as a distinctive morphotype and, 
to indicate the tip-splitting character of the growth, it was denoted T morphotype. In the 
next section we describe in some detail the observations of Ben- Jacob et al. and Matsushita 
et al. . Additional studies of branching colonial growth are reported by Matsuyama et al. 



ggj and Mendelson and Salhi [||]. 

All the various strains reported in the studies quoted above exhibit branching patterns 
during growth on a poor substrate. Drawing on the analogy with diffusive patterning in 



non-living systems ||45|-^, we can state that complex patterns are expected. The cellular 
reproduction rate that determines the growth rate of the colony is limited by the level of 
nutrients available for the cells. The latter is limited by the diffusion of nutrients towards the 
colony (for low nutrient substrate). Hence colony growth under certain conditions should be 
similar to diffusion limited growth in non-living systems as mentioned above PTLESf. The 



study of diffusive patterning in non-living systems teaches us that the diffusion field drives 
the system towards decorated (on many length scales) irregular fractal shapes. Indeed, 
bacterial colonies can develop patterns reminiscent of those observed during growth in non- 
living systems. But, this is certainly not the end of the story. The colonies exhibit a richer 
behavior. This, ultimately, is a reflection of the additional levels of complexity involved 
when the building blocks of the colonies, the bacteria, are themselves living systems. We 
now start to reveal the cybernetic processes (communication, regulation and control) which 
are part of the colonial adaptive self-organization, and their determination of the interaction 
between genetic information and biophysical behavior. 

How should one approach the modeling of the complex bacterial patterning? With 
present computational power it is natural to use computer models as a main tool in the study 
of complex systems. However, one must be careful not to be trapped in the "reminiscence 



syndrome", described by J. D. Cowan as the tendency to devise a set of rules which 



will mimic some aspect of the observed phenomena and then, to quote J. D. Cowan "They 
say: 'Look, isn't this reminiscent of a biological or physical phenomenon!' They jump in 
right away as if it's a decent model for the phenomenon, and usually of course it's just got 
some accidental features that make it look like something." Yet the reminiscence modeling 
approach has some indirect value. True, doing so does not reveal (directly) the biological 
functions and behavior. However, it does reflect understanding of geometrical and temporal 
features of the patterns, which indirectly might help in revealing the underlying biological 
principles. Another extreme is the "realistic modeling" approach, where one constructs an 
algorithm that includes in details all the known biological facts about the system. Such an 
approach sets a trajectory of ever including more and more details (vs. generalized features). 
The model keeps evolving to include so many details that it loses any predictive power. 
Here we try to promote another approach - the "generic modeling" one pO| , p7| , p] , |3T| . We 



seek to elicit, from the experimental observations and the biological knowledge, the generic 
features and basic principles needed to explain the biological behavior and to include these 
features in the model. We will demonstrate that such modeling, with close comparison to 
experimental observations, can be used as a research tool to reveal new understanding of 
the biological systems. 

Generic modeling is not about using sophisticated, as it may, mathematical description 
to dress pre-existing understanding of complex biological behavior. Rather, it means a 
cooperative approach, using existing biological knowledge together with mathematical tools 
and synergetic point of view for complex systems to reach a new understanding (which is 
reflected in the constructed model) of the observed complex phenomena. 

The generic models can yet be grouped into two main categories: 1. Discrete models such 
as the communicating walkers models of Ben- Jacob et al. p7| , |52|j30| and the bions model 



of Kessler and Levine |5yj5^. In this approach, the microorganisms (bacteria in the first 
model and amoebae in second) are represented by discrete, random walking entities (walkers 
and bions, respectively) which can consume nutrients, reproduce, perform random or biased 
movement, and produce or respond to chemicals. The time evolution of the chemicals 



is described by reaction-diffusion equations. 2. Continuous or reaction-diffusion models 



53,|55[. In these models the microorganisms are represented via their 2D density, and a 



react ion- diffusion equation of this density describes their time evolution. This equation is 
coupled to the other react ion- diffusion equations for the chemical fields. In the context of 
branching growth, this idea has been pursued recently by Mimura and Matsushita et al. 



p7|| , Kawasaki et al. and Kitsunezaki |£9|. A summary and critique of this approach 



is provided by Rafols |^ 



Here we describe a new model which includes a lubrication fluid and a model with a 
cutoff, as was proposed by Kessler and Levine [^. We compare the results obtained by 
the various models and the experimental observations. Our main goal is to identify the 
biological and mathematical requirements for branching patterns. We then study the effect 
of nutrients- and signaling- chemotaxis, and conclude that chemotaxis is needed to explain 
the wealth of experimental observations. 



II. EXPERIMENTAL RESULTS 

Several strains of bacteria were reported to produce tip-splitting branched patterns under 
conditions of low level of nutrient. We describe here the experimental results of Ben- Jacob 
et al. p5| , ^^ , |62| - working with colonies of Paenibacillus dendritiformis var. dendron 
(T morphotype) - and Matsushita, Fujikawa, Matsuyama and coworkers |3^j23|j3^j26|j40|l - 
working with colonies of Bacillus subtilis. 



A. Growth patterns of T morphotype 

a. Macroscopic observations All manner of patterns are exhibited by T morphotype 
as the growth conditions are varied. An example of branching pattern is shown in figure 
|I|. This kaleidoscope of shapes may be grouped into a number of "essential" patterns. For 
intermediate agar concentrations (about 1.5% - l.bg in 100m/), at very high peptone levels 
(above IQg/l) the patterns are compact (Fig. ^). At somewhat lower but still high peptone 
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levels (about 5-lOg/l) the patterns, reminiscent of viscous fingering patterns in Hele-Shaw 
devices exhibit quite pronounced radial symmetry and may be characterized as dense 
fingers (Fig. 0b), each finger being much wider than the distance between fingers. For 
intermediate peptone levels, branching patterns with lower fractal dimension (reminiscent 
of electro-chemical deposition) are observed (Fig. ^c). The patterns are "bushy", with 
branch width smaller than the distance between branches. As the peptone level is lowered, 
the patterns become more ramified and fractal-like. Surprisingly, at even lower peptone 
levels (below 0.25g/l for 2% agar concentration) the colonies revert to organized structures: 
fine branches forming a well defined global envelope. We characterize these patterns as fine 
radial branches (Fig. ^). For extremely low peptone levels (below O.lg/l), the colonies 
lose the fine radial structure and again exhibit fractal patterns (Fig. |^a). For high agar 
concentration the branches are very thin (Fig. ^d). 

At high agar concentration and intermediate peptone levels the colonies display a struc- 
ture of concentric rings superimposed on a branching pattern (Fig. ^). At high agar 
concentration and very high peptone levels the colonies display a structure of concentric 
rings in a compact colony (Fig. ^). At high agar concentrations the branches also exhibit 
a global twist with the same handedness, as shown in figure Similar observations dur- 



ing growth of other bacterial strains have been reported by Matsuyama et al. |^ , |26| . We 
referred to such growth patterns as having weak chirality |]52|,p|. 

A closer look at an individual branch (Fig. ^) reveals a phenomenon of density variations 
within the branches. These 3-dimensional structures arise from accumulation of cells in 
layers. The aggregates can form spots and ridges which are either scattered randomly, 
ordered in rows, or organized in a leaf-veins-like structure. The aggregates are not frozen; 
the cells in them are motile and the aggregates are dynamically maintained. 

At the other extreme, of very soft agar (0.5% and below), the T morphotype does not 
exhibit branching patterns. Instead, the growth is compact with density. In the range of 
0.5%-!% agar concentration the colonies typically have a shape of many arms stars. 
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b. Microscopic observations Under the microscope, bacterial cells are seen to perform 
a random- walk- like movement in a layer of fluid on the agar surface. This wetting fluid is 



assumed to be excreted by the cells and/or drawn by the cells from the agar p7| , |62[] . The 
cellular movement is confined to this fluid; isolated cells spotted on the agar surface do not 
move. This is an important observation as we discuss later when formulating the models. 
The fluid's boundary thus defines a local boundary for the branch. Whenever the cells are 
active, the boundary propagates slowly, as a result of the cellular movement and production 
of wetting fluid. 

At very low agar concentrations (below 0.5%) the bacteria swim inside the agar and not 
on its surface. Between 0.5% and 1% agar concentration some of the bacteria move on the 
surface and some inside the agar. 

The observations reveal also that the cells are active at the outer parts of the colony, while 
closer to the center the cells are stationary (do not move) and some of them sporulate (form 
spores). It is known that certain bacteria respond to adverse growth conditions by entering 
a spore stage until more favorable growth conditions return. Such spores are metabolically 
inert and exhibit a marked resistance to the lethal effects if heat, drying, freezing, deleterious 
chemicals, and radiation. 

B. Morphology selection, morphology diagram and velocity-pattern correlations 

The emerging understanding of pattern determination in non-living includes the concepts 
of morphology diagram, morphology selection, morphology velocity correlations and mor- 
phology transitions [|^. In short, the patterns formed in many evolving azoic (non living) 



systems may often be grouped into a small number of "essential shapes" or morphologies 
each representing a dominance of a different underlying effect. If each morphology is ob- 
served over a range of growth conditions, a morphology diagram may exist. The existence 
of a morphology diagram implies the existence of a morphology selection principle and vice 
versa. Ben- Jacob et al. proposed the existence of a new morphology selection principle: the 
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principle of the fastest growing morphology |6^j4^, a principle which should be applicable 



for a wide range of growth conditions. In general, if more than one morphology is a possible 
solution, only the fastest growing one is nonlinearly stable and will be observed, that is, 
selected. 

The new selection principle implies that the average velocity is an appropriate response 
function for describing the growth processes and hence should be correlated with the geo- 
metrical character of the growth. In other words, for each regime (essential shape) in the 
morphology diagram, there is a characteristic functional dependence of the velocity on the 
growth parameters. At the boundaries between the regimes there is either discontinuity in 
the velocity (first order-like transition) or in its slope (second order-like transition). 

At present, there is some evidence for the existence of the new selection principle in 
non-living systems. The new principle might also be valid for pattern determination dur- 
ing colonial development in bacteria pT|,|3|. The bacterial patterns may be grouped into 



a small number of "essential shapes", each observed over a range of growth conditions 
PB| , |5^P3| , |^jr8| . To prove this hypothesis, the next step would be to demonstrate the 
velocity-pattern correlation during colonial growth. 

A plot of the growth velocity as a function of nutrient level for 1.5% agar concentration 



is shown in Ref. [^]. For the presented range of peptone levels it was found that the velocity 
shows three distinct regimes of response, each corresponding also to a distinct morphology 
(the fine radial branches, branching patterns and dense fingers), as was predicted for non- 
living systems. The change in velocity suggests that the switching between morphologies is 
indeed a real morphology transition and not a simple cross-over (see Ref. ||4^). The transition 



at low peptone level (between fine radial branches and branching structure) might be a first 
order morphology transition, i.e. a transition characterized by a jump in the velocity and 
hysteresis. The transition at the higher peptone level (from branches to dense fingering) 
seems to be second-order-like. These observations of velocity-pattern correlations strongly 
support the existence of a morphology selection principle which determines the selected 
colonial morphology for a given morphotype. 
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In non-evolving (equilibrium) systems there is a phenomenon of critical fluctuations when 
the system is kept at the transition point between two phases. At that point the system con- 
sists of a mixture of the two phases. In Ref. it was shown that an analogous phenomenon 
exists in evolving non-living systems and explained that this fact provides additional support 
for the idea of morphology transitions. Figure ^ shows patterns exhibited by colonies grown 
at "critical" peptone levels, where transitions between two morphologies occur. Similarly, 
for the fluctuations displayed by non-living systems, we observe a combination of the mor- 
phologies characterizing the patterns above and below the critical point. These observations 
provide additional support for the relevance of the concepts of morphology selection and 
morphology transition to colonial development. 



C. Growth patterns of B. subtilis 



Matsushita and coworkers p3|,|23|,p4|,p6|,^ studied the colonial branching patterns and 



morphology diagram of the bacteria specie B. subtilis OG-01. A detailed summary of these 



observations is provided by Rafols ||6^. A typical morphology diagram is shown in figure 
p. Note that here the x-axis is the inverse agar concentration and the y-axis is the nutrient 
level. These bacteria are not efficient in producing a lubricating fluid, hence above about 
0.8% agar concentration they can not move on the agar surface: Under such conditions 
and low level of nutrients (below 1 g/1 peptone), DLA-like patterns are observed. As the 
level of nutrients is increased, the patterns become compact, with a cellular structure at the 
interface. 

For low agar concentrations (below 0.5%, so that the bacteria can move) and low level of 
nutrients, dense branching patterns are observed. These patterns are replaced by compact 
growth for higher levels of nutrients. Beautiful patterns of concentric rings imposed on a 
dense branching growth (Fig. 10) are observed at high levels of nutrients and intermediate 



agar concentric (about 0.75%). For more details see Ref pO 



The different morphologies correspond to different growth velocities; DLA-like patterns 
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grow in about a month, compact patterns at intermediate concentrations of agar grow in 
about a week, dense branching patterns and patterns of concentric rings grow in few days, 
and compact patterns at low concentrations of agar grow in half a day. From this we learn 
that indeed the growth velocity of the various morphologies is very different. Moreover, it 
seems that different bacterial movement mechanisms correspond to the different regimes. 
Thus we expect a real transition between the various regimes in the morphology diagram 
rather than a simple cross-over. Therefore, if velocity as function of the growth parameters 
is to be measured, it will probably show a jump in the velocity or its slope (first or second 
morphology transitions, respectively). 



III. BIOLOGICAL BACKGROUND 

Clearly we cannot begin to encompass all the biological background. Thus we will de- 
scribe here, based on our previous experience, only the most relevant information for the 
understanding and modeling of the observed colonial patterning. 



A. Bacterial Surface Translocation 



The most widely studied mechanism used by bacteria for movement is swimming with 
flagella , but other mechanisms exist as well . Most common types of bacterial move- 



ments are categorized to be 

• Swimming - Surface translocation produced through the action of flagella. The cells move 
individually and at random in the same manner as flagellated bacteria move in wet mounts 
(i.e., nearly straight runs separated by brief tumbling). Swimming takes place only in suffi- 
ciently thick surface fluid. Microscope observations reveal no organized flow-field pattern. 

• Swarming - Surface translocation produced through the action of flagella, but unlike swim- 
ming, the movement is continuous and regularly follows the long axis of the cell. The cells 
are predominantly aggregated in bundles during their movement, and microscope observa- 
tions reveal flow-field patterns highly organized in whirls and bands. 
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• Gliding - Surface translocation occurring only in non-flagellated bacteria and only when 
in contact with solid surface. In all other respects, gliding is identical to swarming. 



B. Modeling Bacterial Movement 

As for the movement of T morphotype, based on microscope observations of movement 
and electron microscope observations of flagella we identify the movement as swimming. 
Cells tumble about every ^ 1 — 5 depending on external conditions. The speed of the 
bacterium between tumbling events is very sensitive to conditions such as the liquid viscosity, 
temperature and pH level. Typically, it is of the order of 1-10/im/sec. 



Swimming can be approximated by a random walk with variable step size ||6^. At low 
bacterial densities the random walk can be described by a diffusion equation with a diffusion 
coefficient = v"^ tj- = 10^® — lO^^cm^/sec. Low bacterial densities means that the mean 
free path between bacterial collisions Ic is longer than the tumbling length Ix = vt^, thus 
collisions between the bacteria can be neglected. The mean free path (or collision length) is 

{p~i in 3 dimensions 
a 2 in 2 dimensions 

where p is the 3D bacterial density and a is the 2D density - the projection of p on the 
surface. 

At high densities {Ic < It), the collisions cannot be neglected. In attempt to approximate 
the dynamics in those conditions, one may want to consider the time of straight motion to 
be Ic/v instead of ry. Hence Dh depends on the bacterial density to yield 

{vp~^ in 3D 
(2) 
va 2 in 2D 

This approximation is valid under the assumptions that a collision event is identical to a 
tumbling event (abrupt uncorrelated change in direction of motion), that a tumbling event 
is independent of the collisions, and that the speed between such events is not affected by 
their frequency. 
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The assumption that a colhsion event is hke a tumbhng event poses many problems. 
Even if the bacteria do not activate special response to collision it is unrealistic to assume 
that collisions are elastic, or that the fiagella adopt immediately to the new orientation 
which changes during collisions. Thus it is reasonable to assume strong correlation between 
the cell's orientation before collision and the cell's orientation after collision. In addition, 
the orientation after the collision should be biased according with the average direction of 
motion of the surrounding bacteria, as they carry the liquid with them. The important 
parameter is not the collision length Ic but re-orientation time Tr. The re-orientation time 
is the time it takes a bacterium to loose memory of its initial orientation, i.e. the time span 
on which the final orientation has effectively no correlation with the initial orientation. At 
low densities the re-orientation time is equal to the tumbling time tt- As the density 
rises and the collisions become more frequent, decrease, defines the densities above 
which the constant diffusion coefficient = v"^ tt is not a good approximation. It is quite 
possible that these densities are high enough so as to make the velocity and even the type 
of motion dependent on bacterial density, making relation (0) irrelevant. In any case, high 
cellular densities does mean an effective decrease in the diffusion coefficient related to the 
bacterial movement. 

When swimming in an unstirred liquid, very low cellular densities also effect the move- 
ment. The bacteria secrete various materials into the media and some of them, e.g. enzymes 
and other polymers, change significantly the physical properties of the liquid making it more 
suitable for bacterial swimming. The secretion of these materials depend on cellular density, 
thus at not-too-high densities the speed of swimming rise with the cellular density. Hence 
the diffusion coefficient related to the bacterial movement should be a non-monotonic func- 
tion of the bacterial density. Moreover, the specific functional form might depend on the 
specific bacterial strain. 

In other conditions there is similar but more pronounced effect. On semi-solid surface 
the bacteria cannot swim at all inside the agar and they have to produce their own layer of 
liquid to swim in it. A single T bacterium on the agar surface cannot produce enough fluid 
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to swim in it, thus the bacteria cannot break out of the layer fluid and the branches of a T 
colony can be defined by this fluid. Whenever bacteria enter the shallower parts of the layer, 
at the edge of the branch, they become sluggish, indicating that the depth of the layer effects 
the bacterial movement. It can be argued (see section [IV D| ) that in such cases the bacterial 
speed is related to the bacterial density by a power law (at least in low densities). Not only 
the diffusion coefficient related to the bacterial movement is a non-monotonic function of 
the bacterial density (as in a liquid agar), but it is also vanishes for extremely low densities. 
In this case it is clear that the specific functional form depend on the specific bacterial strain 
{B. subtilis, for example, cannot move at all under such conditions). 



C. Chemotaxis 

Chemotaxis means changes in the movement of the cell in response to a gradient of 
certain chemical fields p7|- |69| , |66| . The movement is biased along the gradient either in the 



gradient direction or in the opposite direction. Bacteria are too short to estimate spatial 
gradients of the chemical by simply comparing concentrations at different locations on their 
membrane (but see Ref. for different view). They deduce the spatial gradients by 



calculating temporal derivatives along their path. It is known that E. coli, for example, can 
compare successive measurements over a time interval of 3 seconds. The actual chemotaxis 
in swimming bacteria is implemented by decreasing the tumbling frequency as cells swim up 
the gradient of the attractant or down the gradient of repellent. Thus the straight runs are 
important for gradient perception and the tumbling timing is important for the response to 
this gradient. 

Usually chemotactic response means a response to an externally produced field like in 
the case of chemotaxis towards food. However, the chemotactic response can be also to 
a field produced directly or indirectly by the bacterial cells. We will refer to this case as 
chemotactic signaling. 

Chemotaxis towards high concentration of nutrients is a well studied phenomenon in 
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bacteria |6^J7T|. When the center of a soft agar plate (0.35% agar concentration) is in- 
oculated with cells capable of chemotaxis, distinct circular bands of bacterial cells become 
visible after a few hours of incubation. In fact, these patterns were used as semi-quantitative 
indicators of chemotactic response ||7I|. Genetic experiments showed that the creation of 
each of those bands depends solely on the chemotactic response to a single chemical in the 
substrate (these chemicals are usually metabolizable, but even cells that have lost the abil- 
ity to metabolize a certain chemical form bands, as long as they are attracted to it [|67[] ). 
Berg et al. showed that the bacteria realize chemotactic response by modulating the 



periods between tumbling events - they decrease the probability of tumbling when moving 
in a preferred direction along the chemical gradient. This makes a bias in the random walk 
which result in a mean drift of the bacteria in the desired direction, a drift that can be as 
large as f/10. 

Bacteria sense the local concentration C of a chemical via membrane receptors binding 
the chemical's molecules |6^J6^. The cell measures the concentration by calculating the 



relative number of occupied receptors -j^^^,, where No and Nf are the number of occupied 
and free receptors respectively. For a given chemical C, No is determined by two character- 
istic times: the mean time of a receptor occupation - Tq, a constant determined by internal 
cellular processes - and the mean time lapse when the receptor is free {tj). Since Tf is 
inversely proportional to the concentration of the chemical (with the proportion coefficient 
determined by the receptor affinity to the chemical), we get: 

No ^ To ^ C 

Nf + No Tf + To K+C ' ^ ^ 

Where K = {C Tf) /to is constant. It is crucial to note that when estimating gradients of 
chemicals, the cells actually measure changes in the receptors' occupancy No/ {No + Nf) and 
not in the concentration itself. Using Eq. (|^) we obtain: 

d_ ( No \ ^ K dC 

dx \No + Nf) (i^ + C)2 dx' ^ ' 

This means that the chemical gradient times a factor K/ {K + C)^ is measured. This de- 
pendence in known as the "receptor law" For very high concentration the chemotaxis 
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response vanishes due to saturation of the receptors. The chemotactic response also vanishes 
at the opposite hmit of small concentration, as the concentration reception is masked by 
external and internal noises. This effect is not included in the receptor law, which should 
be changed accordingly. 

The receptor law is needed to explain the the bands reported by Adler [^jO • It can be 
shown that linear chemotactic response to a nutrient cannot produce such bands. A non- 
linear response like the "receptor law" must be included for the bands to form. Moreover, 
high concentration of the attractant represses both the strength of the chemotactic response 



and the velocity of the expanding band [74|. These observations are accounted for by 



the "receptor law" for chemotactic response if one assumes that the average gradient sensed 



by the cells is proportional to the initial concentration of the chemical ||67| , [74 
The bacterial flux due to chemotaxis can be described by 

Jchem = C(^)X(C)VC (5) 

where x(C)VC is the gradient sensed by the cell (with x(C) having the units of 1 over 
chemical's concentration) and is the bacterial response to the sensed gradient (having 
the same units as a diffusion coefficient). xiC) is usually taken to be either constant or the 
"receptor law". 

The function of the bacterial response C, is positive for attractive chemotaxis (movement 
towards high concentrations) and negative for repulsive chemotaxis. If the movement is 
in liquid and at low bacterial densities, then \C{(^)\ oc aDf,. In a lubrication fluid which 
effect the bacterial movement, the chemotaxis is effected in the same way the diffusion is; 
\C{a,l)\^aD,{a,l). 

In the case of high bacterial densities, collisions between bacteria can disrupt both the 
perception of chemical gradient and the bacterial response. As the collisions prevent the 
bacteria from moving on a straight line between tumbling events, the effective response to 
chemotaxis is reduced. 
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D. Food Consumption, Reproduction and Starvation 



The T bacteria, like most bacteria, reproduce by fission of tlie cell into two daughter 
cells which are practically identical to the mother cell. The crucial step in the cell division 
is the replication of the genetic material and its sharing between the daughter cells. Haste 
replication of DNA might lead to many errors - most organisms limit the rate of replication to 
about 1000 bases per second. Thus the reproduction must take at least minimal reproduction 
time T/j. This reproduction time tr is about 25mm in Bacilli. 

For reproduction, as well as for movement and other metabolic processes, bacteria and all 
other organisms need influx of energy. Any organism which does not get its energy directly 
from sunlight (by photo-synthesis) needs an external supply of food. In the patterning 
experiments the bacteria eat nutrient from the agar. As long as there is enough nutrient 
and no significant amount of toxic materials, food is consumed (for cell replication and 
internal processes) at maximal rate fi^- To estimate fic we assume that a bacterium needs 
to consume an amount of food Cr of about 3 x IQ^^'^g. It is 3 times its weight - one quanta 
for doubling body mass, one quanta used for movement and all other metabolic processes 
during the reproduction time tr, and one quanta is for the reduced entropy of making 
organized cell out of food. Hence VL^. is about Ijgj sec {Ifg = 10~^^gram). 

If nutrient is deficient for a long enough period of time, the T cells may enter a special 
stationary state - a state of a spore - which enables them to survive much longer without 
food. The bacterial cells employ very complex mechanisms tailored for the process of sporu- 
lation. They stop normal activity - like movement - and use all their internal reserves to 
metamorphose from an active volatile cell to a sedentary durable 'seed'. While the spores 
themselves do not emit any chemicals (as they have no metabolism) , the pre-spores (sporu- 
lating cells, see Fig. ^) do not move and emit a very wide range of waste materials, some of 
which unique to the sporulating cell. These emitted chemicals might be used by other cells 
as a signal carrying information about the conditions at the location of the pre-spores. Ben- 



Jacob et al. p7| , |62| , [75| suggested that such materials are repelling the bacteria ('repulsive 
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chemotactic signaling') as if they escape a dangerous location. 

When bacteria are grown in a petri dish, nutrients are usually provided by adding pep- 
tone, a mixture including all the amino acids and sugars as source of carbon. Bacteria which 
are not defective in synthesis of any amino acid can grow also on a minimal agar in which 
a single source of carbon and no amino acids are provided. Such growth might seem to 
be easier to model as the growth is limited by the diffusion of a single chemical. However, 
during growth on minimal agar there is usually a higher rate of waste products accumu- 
lation, introducing other complications into the model. Moreover many of our strains are 
auxotrophic i.e. defective in synthesis of some amino acids and need an external supply of 
it. Providing the bacteria with these amino acids and only a single carbon source might 
pose us the question as to what is the limiting factor in the growth of the bacteria. For all 
those reasons we prefer to use peptone as nutrient source. 

We said that if there is ample supply of food, bacteria reproduce in a maximal rate of one 
division in Tpt- If the available amount of food is limited, bacteria consume the maximum 
amount of food they can. In the limit of low bacterial density, the available amount of food 
over the tumbling time tt is the food contained in the area TT\fD\J\.., where Df, and D„ are 
the diffusion coefficients of the bacteria and the food, respectively. Hence the rate of food 



consumption is given by n^D^yDn (weather D^, is constant or not). 

In a continuous model, reproduction of bacteria translate to a growth term of the bacterial 
density which is a times the eating rate per bacteria. In the limit of high nutrient it is ct/tr, 
and in the limit of low nutrient it is proportional to no. This brings to mind Michaelis- 
Menten law [^] of j^^na with K, 7 constants. Many authors take only the low nutrient 
limit of this expression, Kna, although it is not biologically established that the bacteria in 
the experiments are limited by the availability of food and not by their maximal consumption 
rate. 
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IV. REACTION-DIFFUSION MODELS 



In this section we deal with continuous, reaction-diffusion models for bacterial growth. 



The models under study are due to Fisher and Kolmogorov |76.77|, Kessler and Levine [61 
Kitsunezaki Kawasaki et al. and Mimura et al. The models are all two- 

dimensional (2D), with 6(x, t) denoting the density of bacteria projected on a 2D plane and 
n(x, t) is the 2D nutrient density. The equations for the various models will be written in 
dimensionless units, and the reader is referred to the Appendix for a discussion about the 
relations with real units. 



In general, the rate of change of the bacteria density can be described by i73|| : 

dh 

— = movement + ''birth" — ''death" (6) 
at 

As discussed in Sec. JTT|, the movement of bacteria consists of various possible mechanisms, of 
which we will concentrate on swimming, so that the motion is described as diffusion (either 
linear or non-linear). The "birth" term in Eqn. ^ corresponds to bacterial reproduction, 
which depends on the supply of nutrients. The "death" term represents the transition of 
bacteria into a non moving state. 

A. The Fisher-Kolmogorove equation 

Mathematically, the above description is usually written as a reaction diffusion equation, 
for which the canonic example is the Fisher-Kolmogorove equation 0,|7^] (without a death 
term): 

= D.V'b + b{l - b) (7) 

This equation was originally presented to describe the spread of mutants in a population. 
We will use it here as a starting point for our discussion of colonial development. In this 
context, Df) is the diffusion coefficient describing the bacterial movement, and the reaction 
term / (6) = 6(1 — 6)) describes both the growth and "death" of bacteria. The function f{b) 



is depicted in Fig. |[ Eqn. ^ has two homogeneous solutions, a stable solution 6 = 1 and 
an unstable solution 6 = 0. These solutions correspond to the two extrema of the potential 
$ = — / f{b)db, which can be thought of as a Landau-Ginzburg free energy density (see 
Fig. Thus, we can study the propagation of the stable state (inside the colony) into the 
unstable one (outside the colony). It can be shown |7^-§3| that in ID there is a unique 
selected velocity of the front, chosen according to the properties of the fixed point far ahead 
at infinity. 

In two dimensions, the propagation is in the form of a compact (as opposed to branching) 
growth. That is, there is no MuUins-Sekerka instability [p| , ^^ -|86i- Iii the case of such 



instability, a small bump in a flat interface will have a higher velocity than the rest of the 
front and will therefore over-grow. Here, however, this will not happen, because the front 
velocity is determined by the fixed point at infinity, rather than by local properties of the 
front. 

The Fisher-Kolmogorove model is an appropriate description of the growth when bacteria 
are grown under nutrient-rich conditions, in which case the growth dynamics is not limited 
by the supply of food P^. Here we are interested in the opposite case, where nutrient 



supply is limited. For a more realistic description of colonial development on a nutrient- 
poor surface, we must take into account the interaction of bacteria with the nutrient field 
n{'x,t). In the simplest case, this is described by the Diffusive Fisher-Kolmogorove equation 

8h 

^^=D,V% + f{b,n) (8) 
^ = D^V'n-r^f{b,n) (9) 

where 77 > is the conversion ratio of food into bacteria (3 picogram per bacteria, see Sec. 
pnp . In this model the food consumed by the bacteria is reduced from the food field and 
converted into bacteria. The shape of a one- dimensional front obtained for f{b, n) = bn is 
depicted in Fig. [1^. It turns out that this model (as the original Fisher-Kolmogorove model) 
has a selected front velocity, determined by the conditions at infinity, and therefore does 
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not exhibit diffusive (Mullins-Sekerka) instability, and a two dimensional growth is compact 
rather than branched. An additional way of understanding why compact growth is obtained 
even under poor nutrient levels is to note that the model, as it is, does not impose a minimal 
bacteria density in the colony, so that the density can be adjusted according to the initial 
food level. In the real biological system, however, some minimal density is required in order 
to create the lubricating fluid (needed for movement), and so compact growth usually is not 
possible. 

For further understanding the requirements for a branching pattern, let us recall the case 
of solidification from an under-cooled melt, which exhibits branching patterns ||47|j48|l . Our 
bacterial density b corresponds to the order parameter in the phase-field model description 
of solidification ||88| , |89[| , whereas the diffusion of food is analogous to the diffusion of heat 
away from the solid-liquid interface. In the case of solidification, the Landau-Ginzburg free 
energy is a tilted double- well (see Fig. [Tl|). The meta-stable state corresponds to the liquid 
phase. In the diffusive Fisher-Kolmogorove case, the analog of the liquid state, i.e. the 6 = 
state, is unstable. Thus, according to the solidification case, if we modify the model and 
turn the 6 = state into a meta-stable one, this can lead to branching growth. 



B. A cutoff in the reaction term 

In the case of bacteria, there is a feature of the system that might have a similar effect to 
the meta-stability in solidification. This is the discreteness of bacteria, for which the contin- 
uous description is not always valid. Kessler and Levine |^ argue that when describing a 



discrete system using continuous models, a cutoff near the fixed point must be imposed, i.e. 
the reaction term must be set to zero when the (bacterial-) density is below some threshold. 
They have shown that inclusion of such a cutoff leads to a Mullins-Sekerka instability, and 
branching patterns may appear when a death term is also included, as explained below. 
To show this, we consider the diffusive Fisher equations with a cutoff: 

— = DbW% + bnQ{b - /3) (10) 
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^ = - nbQ{b - (3) (11) 

where f5 is the threshold density for growth, and G is the Heaviside step function. The 
food consumption term is of the form f{n,b) = nb, which is the widely used low- nutrient 
approximation for the Michaelis-Menten law (Sec. |1II D| ). The value of the threshold for 
bacterial growth, f3, is taken to be one bacterium per 1-10 /im^. Note that this corresponds 
to the case where the distance between bacteria is of the order of the length between tumbles. 
Fig. O depicts the reaction term and the Landau-Ginzburg free energy for this model. As 



can be seen in Fig. |T3], an instability of the front indeed appears, and the compact growth 
pattern has a surface broken by "fjords". However, the Mullins-Sekerka instability is not 
sufficient to produce branches. The emerging dips soon "heal", so that branches are not 
formed. One way to obtain branching growth is to add a "death" term to the model, thus: 

db 

— = DhW% + bnQ{b - /3) - nb (12) 

^ = V'n - nbQ{b - (3) (13) 

1^.6 (14) 

where /i is the rate of bacterial differentiation into non-moving state, and s(x, t) is the 
density of "frozen" bacteria. This modified model exhibits distinct branching patterns, as 



seen on Fig. |T^. The explanation for this effect lies in the fact that now, with a death term 
present, bacteria left behind the propagating front become non-motile ("dead"). They are 
unable to move and close the "fjords" , thus allowing real branches to form. 



C. Reaction-Diffusion with Lubrication 

We have so far ignored the effect of the lubricating field on the motion of the bacteria. 
We present here a new model which incorporates an additional field that describes the 
lubricating fiuid. The field, denoted as /, is the local height of the lubrication fiuid on 
the agar surface. Its dynamics is governed by a reaction diffusion equation. The bacterial 
diffusion is coupled to this field. 
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The dynamics of the field is given by: 

dl 

- = -VJi + rbn{lma. -I) -XI (15) 

where Ji is the fiuid fiux (to be discussed) , F is the production rate and A is the absorption 
rate of the fluid by the agar. 

We assume that the fluid production depends on the bacterial density. As the production 
of lubrication probably demands substantial energy, it should also depend on the nutrients 
level. We assume that the absorption of fluid into the agar depends on the local amount 
of fluid (i.e. the height of the fluid layer). In this model the production depends linearly 
on the concentrations of both the bacteria and the nutrients. The production term cannot 
become negative as the lubrication height cannot exceed Imax- 

The lubrication fluid flows by diffusion and by convection caused by bacterial motion. 
A simple description of the convection is that each bacterium drags along its movement the 
fluid surrounding it. 

Ji^-Dii^m+jJb (16) 

— * 

where Di is the lubrication diffusion coefficient, Jf, is the bacterial ffux and j is the amount 
of fluid dragged by each bacterium. The diffusion term of the fluid depends on the height 
of the fluid to the power 77. The nonlinearity causes the fluid to have a sharp boundary at 
its front, as is observed in the experiments of bacterial colonies development. 

We now turn to the effect of the lubrication on bacterial diffusion. An increase in the 
amount of lubrication decreases the friction between the bacteria and the agar surface. The 
term 'friction' is used here in a very loose manner to represent the total effect of any force 
or process that slows down the bacteria. It might include, for example, the drag which acts 
on a body moving in shallow layer of viscous fluid. It might include the probability that a 
flagellum will adhere or get tangled with the polymers of the agar. As the bacterial motion 
is over-damped, the local speed of the bacteria (or the mean step length, when neglecting 
collisions between bacteria) is proportional to the self-generated propulsion force divided 
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by the friction. It can be shown that variation of the step length leads to variation of 
the diffusion coefficient, with the diffusion coefficient proportional to the step length to the 
power of two. We assume that the friction is inversely related to the local lubrication height 
through some power law: friction ~ P and 7 < 0. Following our arguments the bacterial 
flux is: 

Jb = -Dbl-^^Vb (17) 

For the complete model we took simple bacterial growth and death terms. The model is: 
dh 

— = DbV{l-^^Vb)+bn- i^b (18) 



V^n — bn 



dn 

'dt 
81 

- = V{DirVl + jDJ-^^Vb) + TbniLa. - I) - M 
ds , 



(19) 



For the initial condition, we set: 



n{x,t) = no (20) 
6(x,t) = 6o(x) 

where hq is the initial concentration of nutrients and 60 (x) is the initial bacterial distribution. 
In the simulations, &o(x) is localized at the center. 

Preliminary results show that the model can reproduce branching patterns (Figure |15|) . 



At low values of the absorption rate the model exhibits dense fingers. At higher values of 
the absorption rate the model exhibits finer branches. We obtain finer branches also if we 
change other parameters that effectively decrease the amount of lubrication. We can relate 
these conditions to high agar concentration. In this model, as in the non-linear diffusion 
model described below, the bacterial field and the lubricating field have a front "wall" with 
compact support (Figure 0). 
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D. Non-linear diffusion 



It is possible to introduce a simplified model, where the fluid field is not included, and 
is replaced by a density- dependent diffusion coefficient for the bacteria Db ~ b'' |]90|,|9T 



For this purpose, a few assumptions are needed about the dynamics at low bacterial and 
lubrication density: 

• The production of lubricant is proportional to the bacterial density to the power a > 0. 

• There is a sink in the equation for the time evolution of the lubrication field, e.g. 
absorption of the lubricant into the agar. This sink is proportional to the lubrication 
density to the power /3 > 0. 

• Over the bacterial length scale, the two processes above are much faster than the 
diffusion process, so the lubrication density is proportional to the bacterial density to 
the power of (3 /a. 

• The friction is proportional to the lubrication density to the power 7 < 0. 

Given the above assumptions, the lubrication field can be removed from the dynamics and 
be replaced by a density dependent diffusion coefficient. This coefficient is proportional to 
the bacterial density to the power k = — 27/^/0; > 
A model of this type is offered by Kitsunezaki 

db 

= \/(Dob''Vb)+nb- fib (21) 

ot 

^ = V^n - bn (22) 
dt ^ ' 

(23) 

For A; > the ID model gives rise to a front "wall", with compact support (i.e. 6 = 
outside a finite domain, see Fig. |T^). For k > \ this wall has an infinite slope (Fig. p!8[ ). 



The propagation velocity in this case is determined by the condition at the front, not at 
infinity |9^j9^. We therefore expect a MuUins-Sekerka instability in 2D (as is claimed in 
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). Indeed, the model exhibits branching patterns for suitable parameter values and 
initial conditions (Fig. |1^). Note, that the compact support exhibited by this model (that 
is, the abrupt disappearance of bacterial presence outside the colony boundary) is much 
more in accordance with experimental observations than the long "tail" of bacterial density 
appearing in the Fisher-Kolmogorove case (recall Fig. |TU|). 

Another state-dependent diffusion coefficient was proposed and studied by Kawasaki 
et al. which took ~ nh. They justify this form by the observation that bacteria 



are active mostly at the edge of the colony - the only area where there is high amount of 
bacteria and food. Their model, too, exhibits branching shapes (see Fig. ^Oj). This is due 
to the h dependence of the diffusion coefficient, which leads to front instability, just as in 
the Kitsunezaki model. The fact that -Dj, also depends on n prevents bacteria inside the 
colony from moving - and closing the dips created by the instability. In this way, branches 
are created without a need for a death term. A similar mechanism of "food" dependent 
diffusion coefficient was used by Tu et al. [^, who describe a mean-field model for DLA. 
Their model, too, does not include a death term yet produces branching patterns. 

E. A meta-stable reaction term 

As mentioned earlier, meta-stability of the growth term can lead to branching patterns. 



Mimura et al. ||56| , |57| have studied the following model, for which 6 = is a meta-stable 
fixed point: 

^ = D,V'6 + ehn - ^ ^ (24) 

dt ^ (6 + l)(n + l) ^ ' 

^ = y^ri- nh (25) 
dt ^ ' 

ds ^ fib 

dt (6+l)(n + l) ^ 

(see the Appendix for relations to real units). In this model, Db does not depend on bacterial 
density, and it's value is said to vary with agar concentration. The key feature of this model 
is that the total bacteria growth term (multiplication minus inactivation) , depicted in Fig. 
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2T| , gives a meta-stable state at b=0. This means that in order to initiate bacterial growth, 
a threshold value of b* = — 1 must be reached (this value corresponds to bacterial 

density of the order of 0.1 bacteria/yum^). This feature contradicts the observations, that 
even a small number of bacteria, inoculated on a substrate, will multiply and later begin to 
move. This does not imply that the model is incorrect. A possible interpretation is that 
the field b actually describes a combination of lubricant + bacteria. In this case, however, 
as we have explained before, we would expect the model to exhibit a non-linear diffusion. 
Hence we believe that this model might provide a better description of the bacteria if the 
diffusion term is replaced with a non-linear diffusion term. The model as it is exhibits 
various branching patterns, patterns of concentric rings and compact growth (Fig. ^31) . 

Mimura et al. argue that the model captures the experimental morphology diagram they 
observed. This is a very crucial point. If indeed the above claim is correct, it implies that 
the observed patterns can be reproduced with no need for additional biological features. 



However, using the discrete communicating walkers model Ben- Jacob et al. have con- 
cluded hat the additional features of chemotactic response have to be included. So, in order 
to check this point we have performed more detailed comparison between the Mimura et al. 
model and the experimental observations. 

First, we consider the DLA-like growth. In this case, the bacteria do not move on the 
agar surface, and the growth is indeed very similar to the DLA algorithm, as was proposed by 
Matsushita et al. ||37| , |33|1 . It is now understood that in a mean- field DLA model the particle 



density (density of bacteria in the present case) can not be described by a diffusion term. 



Instead, it has to be described by a diffusion multiplied by the nutrients field [Q, which 
differs from the linear diffusion in the model of Mimura et al. . Indeed, close inspection of 
the fractal pattern created by the model reveals that it differs from the observed DLA-like 
patterns. 

Another test of the model is the predicted pattern of concentric rings. It has already 
been pointed out by Rafols pD| that the model's pattern differs from the observed one. In 



the experiment, branching growth slows down. The branches become wider and growth 
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stops. Then, after bacterial differentiation, a new cycle of branches growth starts with thin 
branches emitted from the stationary wide branches [^. This description differs from the 
model patterning shown in Fig. |2^. 

In Fig. ^ we exhibit results of numerical simulations for various levels of peptone and 
for agar concentration for which concentric rings are observed. The sequence of patterns 
from DLA-like at low peptone to concentric rings at high levels of peptone differ from the 
similar sequence of observed patterns presented in Ref. [|60[] . 

We have also tested the change in patterns as we vary the agar concentration (see Fig. 



251) . When we plot the growth velocity as a function of the agar concentration, it does not 
show a jump in the velocity or its slope. In other words, the model seems to exhibit a simple 
crossover between the patterns rather then a morphology transition as the observations seem 
to indicate. 

The above results lead us to conclude that the Mimura et al. model does capture some 
of the observed branching patterns, yet the complete description of the observations requires 
additional features to be included. Specifically, we propose to include nonlinear diffusion. 
We also believe that chemotactic response does play an important role for poor growth 
conditions. To further test this conclusion, we present in the next section a study of patterns 
produced by the reaction-diffusion models when chemotactic responses are included. 

V. INCORPORATION OF CHEMOTACTIC SIGNALING IN THE 

CONTINUOUS MODELS 

So far we have seen several models for the branching colonies, each with its own mecha- 
nism of diffusive instability, which produces patterns resembling the observed ones. Is there 
a way to distinguish between the models so as to find out what are the biological features 
underlying the diffusive instability? We rely on results from an atomistic model - the Com- 
municating Walkers model p7|j9^J75[] - for an indication of what are the biological features 
relevant to the different morphologies. 
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A. Chemotactic-based branching growth: Insights from the Communicating Walkers 

model 



Ben- Jacob et al. argued that for the colonial adaptive self-organization the T morpho- 
type employs several kinds of chemotactic responses. Usually chemotactic response means 
a bias of the movement in response to a gradient of an externally produced field like in the 
case of food chemotaxis. However, it could also be a response to a field produced directly 
or indirectly by the bacterial cells - chemotactic signaling. 

Three kinds of chemotactic responses are suggested to be employed by the T morphotype, 
each dominant in different regime of the morphology diagram. One response is the food 
chemotaxis we have mentioned earlier. According to the "receptor law", it is expected 
to be dominant for some range of nutrient levels (the corresponding levels of peptone are 
determined by the constant K) . The two other kinds of chemotactic responses are signaling 
chemotaxis. One is long-range repulsive chemotaxis where the chemical is secreted by starved 
bacteria at the inner parts of the colony. The second signal is a short-range attractive 
chemotaxis where the chemical is secreted by bacteria at the colony's front, bacteria which 
are immersed in toxic waste products. The length scale of each signal is determined by the 
diffusion constant of the chemical agent and the rate of its spontaneous decomposition. 

Amplification of diffusive Instability Due to Nutrients Chemotaxis: In non-living systems, 
the more ramified patterns (lower fractal dimension) are observed for lower growth velocity. 
Based on growth velocity as function of peptone level and based on growth dynamics, Ben- 
Jacob et al. concluded that in the case of bacterial colonies there is a need for mechanism 
that can both increase the growth velocity and maintain, or even decrease, the fractal 
dimension. Ben- Jacob et al. suggested food chemotaxis to be the required mechanism. It 
provides an outward drift to the cellular movements; thus, it should increase the rate of 
envelope propagation. At the same time, being a response to an external field it should also 
amplify the basic diffusion instability of the nutrients field. Hence, it can support faster 
growth velocity together with a ramified pattern of low fractal dimension. 
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The Communicating Walkers model was used to test the above hypothesis. The sim- 
ulations showed that as expected, the inclusion of food chemotaxis in the model led to a 
considerable increase of the growth velocity without significant change in the fractal dimen- 
sion of the pattern. 

Repulsive Chemotactic Signaling: We focus now on the formation of the fine radial 
branching patterns at low peptone levels. From the study of non-living systems it is known 
that, in the same manner that an external diffusion field leads to the diffusion instability, an 
internal diffusion field will stabilize the growth. It is natural to assume that such a field is 
produced by some sort of chemotactic agent. To regulate the organization of the branches, 
it must be a long-range signal. To result in radial branches it must be a repulsive chemical 
produced by cells at the inner parts of the colony. The must probable candidates are the 
bacteria entering a pre-spore stage due to depletion of nutrient. This proposal has also been 
verified by simulations of the Communicating Walkers model. In the presence of repulsive 
chemotaxis the patterns become much denser with a smooth circular envelope, while the 
branches are thinner and radially oriented. 



B. Results for the continuous models 

We incorporate the effect of chemotaxis in the continuous models by introducing a chemo- 
tactic flux Jchem, which is written (for the case of a chemorep client and a linear diffusion) 



as 



Jchem = -hx{r)Vr (27) 

where r(x, t) is the concentration of chemorepellent and %('") is the chemotactic sensitivity 
to the repellent. In the case of a chemoattractant, e.g. a nutrient, the expression for the 
flux will have an opposite sign. 

Recall that in the case of the "receptor law", the sensitivity x(r) takes the form: 



Xir) - ^ (28) 
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Thus, we obtain the reaction diffusion chemotaxis equation: 

db 

— = _v(-Z}V6-6x(r)Vr) + /(6) (29) 

In addition, one has to write an equation describing the diffusion and the production 
and decomposition of the chemorepellent. This is written as follows: 

f = D^VV + r^s - rirbr - Kr (30) 

where Dr is the diffusion coefficient of the chemorepellent, is the emission rate of repellent 
by pre-spores, fir is the decomposition rate of the repellent by active bacteria, and is the 
rate of self decomposition of the repellent. 

We have tested the effect of food chemotaxis and repulsive chemotaxis in several of 
the aforementioned reaction-diffusion models, and present here results for the Kitsunezaki 
nonlinear diffusion model and the Mimura et al. meta-stable model. 

When treating the Non-linear diffusion model, we must modify the expression for the 
chemotactic flux, in order to incorporate the density dependence of the bacterial movement. 
Thus, similarly to the diffusion coefficient, which is written as = DqB'', we modify the 
chemotactic flux to become: 

Jchen. = C{b)x{r)Vr (31) 

where ({b) = b ■ = b^^^ is the bacterial response to the chemotactic agent. For the linear 
diffusion case, C,{b) degenerates again to b. 

Fig. ^ depicts a pattern developed by the Kitsunezaki model when food chemotaxis is 
included. All of the parameters are the same as in Fig. |T9| (no chemotaxis). Although the 
patterns are very similar, the growth velocity when food chemotaxis is included is about 
twice the velocity in the absence of chemotactic response. In other words, the velocity is 
doubled with no change in the fractal dimension. 



The effect of repulsive chemotactic signaling is demonstrated in Fig. |28| - again with 
otherwise identical parameters to those in Fig. |r^. It can be seen that the previously 
fractal-like shape has turned into a radial branching pattern with a circular envelope. 
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Thus, for the two types of chemotactic response - food and repulsive, the results we 
observe are similar to those obtained for the communicating walkers model. This agreement 
is not surprising, as both models capture the important feature of a lubrication fluid. Recall 
also that the nonlinear diffusion model is the only one which exhibits a sharp front in ID. 

The effect of chemotactic response in the Mimura et al. model is presented in Figs. 



and The chemotactic response was added to a previously DLA-like colony (Fig. ^ a). 
The addition of food chemotaxis turns the colony into a densely branched one, with branches 
much thicker than before. The repulsive chemotaxis makes the branches radially oriented, 
but they become thicker than before. Thus, the effect of chemotactic responses in this model 
differs from the one obtained for both the walker model and the Kitsunezaki model. We 
believe this to stem from the fact that the model does not include nonlinear diffusion to 
capture the effect of lubrication. 



VI. DISCUSSION 

We first briefly reviewed experimental observations of branching patterns in various bac- 
terial strains, under a range of growth conditions. Both the colonial patterns and the optical 
microscope observations of the bacterial dynamics were presented. We have also included 
a brief summary of the known key biological features required, as we think, for successful 
modeling of the growth. 

Our goal in this manuscript was to test reaction-diffusion models. To this end we surveyed 
the reaction-diffusion models for branching growth that we are familiar with. Mathematical 
analysis reveals that a number of different features can lead to instability of a propagating 
front: A density-dependent (i.e. nonlinear) diffusion coefficient; a lower cutoff in the growth 
term, and a meta-stable growth term. For this instability to create pronounced branches, 
a death term must be added to the growth. Such a term prevents bacteria inside the 
colony from moving and "healing" the dips on the surface. Making the diffusion term 
nutrient dependent can lead to a similar effect. The experimental observations provide a 
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clear indication that the bacteria turn immobile, so that a "death" term indeed has to be 
included in the models. 

The fact that different mathematical features can lead to similar (to the eye) branching 
patterns emphasizes how cautious we have to be in modeling the colonial development. True, 
it might be that different bacterial strains develop branching patterns by the employment 
of different biological features (each might correspond to one of the mathematical features) . 
Yet, when we consider a specific bacterial strain, comparing the experimentally observed 
and model's patterns is not sufficient to tell us if indeed the right biological features are 
included in the model. 

What, then, is the right way to tackle the problem? The integrative way: One must 
combine the mathematical knowledge (in this case, what mechanisms lead to instability 
and branching) with an attempt to model the (assumed) generic biological features, and 
comprehensive comparison with experimental results. As an example, the (mathematical) 
mechanism of non- linear diffusion leads to a branching colony, with compact support. We 
must then ask ourselves whether the movement of bacteria could be described in such a way. 
The answer is positive. As we have explained earlier, the lubrication fluid can be modeled 
by a nonlinear diffusion of the bacteria. 

For more critique test of the models, additional aspects of the growth (such as functional 
dependence of the colonial growth velocity on growth conditions, branches size and width 
distributions etc.) have to be compared with the model predictions. One should also compare 
the theory with more involved experimental tests, such as the effect of imposed anisotropy, 
competition between neighboring colonies, and expression of mutants (emergence of sectors) 
in expanding colonies. 

Our conclusion from the study of bacterial branching growth is that the minimal features 
of diffusion, food consumption, reproduction and inactivation are not sufficient to explain 
the complete picture of the observed phenomena. We believe that additional mechanisms 
must be introduced, and propose chemotactic signaling as plausible one. 

This work has dealt with continuous models. Such models are preferable to discrete ones. 
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Each has its advantages and disadvantages. The discrete walkers model, for example, enables 
us to include the valuable feature of internal degrees of freedom, but is computationally 
limited in the number of walkers that can be simulated, and thus its scaling to the real 
problem is somewhat problematic. The best strategy is to employ in parallel both the 
react ion- diffusion and the walkers approach. 

To conclude, despite the difficulties and possible pitfalls, we hope to attract researchers 
to this emerging new endeavor. After all, a significant progress has been made towards 
working out the cybernetic processes (communication, regulation and control) during colo- 
nial development. Yet, many challenging puzzles are waiting to be solved and tantalizing 
phenomena are waiting to be discovered. 
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APPENDIX: FROM REAL UNITS TO DIMENSIONLESS MODELS 

Our goal is to relate the dimensionless equations with biophysical values of the parame- 
ters. The procedure is to set the dimensional units to be the natural scales. 

A. Deriving the equations 

The diffusive Fisher-Kolmogorove equations in dimensional form (section [IV A|) are: 
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db 

— = DtW% + E,nb (32) 

— = DnV^n - Enub 
at 

b and n are the two dimensional concentration of bacteria and of nutrition respectively. 
The field b is measured in units of number of bacteria per area in cm?. The field n is 
measured in units of grams per area in cth? . Experimentally nutrients are usually measured 
in gram /liter. Note that Ig/l corresponds to O.Smg/cm"^. Dh and D„ are the corresponding 
diffusion coefficients in units of am? / sec. Eh is the bacterial reproduction rate in units of 
sec~^ per nutrition concentration. En is the nutrition consumption rate in units of sec~^ 
per bacteria concentration. 

We change the variables to be dimensionless. 

t ^ tT X — >■ xX 

(33) 

b ^ bB n ^ nN 

the new variables are dimensionless and the capitals are their corresponding units. We let 
the temporal and spatial units be the natural scales. The (microscopic) time scale of the 
model is the bacterial reproduction time r^j, so T = r^j. The length scale is the diffusion 
length of the nutrition during reproduction time. The nutrition available for bacteria during 
reproduction time is proportional to the square of the diffusion length. The length unit is: 



X = JDnTR (34) 



After inserting the dimensionless variables into Eq. (|3^) we obtain: 

db 

di 
dn 



Dh/DnV%+{TN)Ehbn (35) 



^ = V^n - {TB)EM 

We define the relative diffusion coefficient D = Db/Dn and impose the following relations: 

TBEn = 1 (36) 
TNEh = 1 (37) 
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We obtain the following dimensionless equations: 

db 

di 
dn 



DV% + hn (38) 

at 



^ = V^n — hn 
dt 



B. Evaluation of the parameters 

We will estimate the values of the parameters Eb, En and Z)„, and derive from them the 
dimensional units. A review of some of the following biological arguments appears in section 

• The bacterial reproduction time, when bacteria grow under optimal conditions, is 
about tr = 25mm. Colonies which exhibit branching patterns grow under limited 
nutrition supply. Therefore the reproduction time will be longer, but in the same 
order of magnitude. We set the time unit to be: 

T = tr = 25mm (39) 

• A typical value for the diffusion coefficient of chemicals in agar is 10~^cm^/sec. So we 
assume that Dn, the diffusion coefficient of the nutrition in the agar, is similar. We 
can find the length unit using 



X = 0.01cm = 100/xm (40) 
The nutrition concentration in the experiments conducted by Ben- Jacob et al. 0, 



Rafols |60] and others was 0.1 — bmg/am?. We set to have a similar value: 

N = Img/cm'^ = 10 x lO'^^g/fi^ = lO'^g/X^ (41) 



The reproduction rate per bacterium (Eq. |32|) is EbNn, where n is the dimensionless 
concentration. The rate is the inverse of the reproduction time, which depends on 
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the nutrition concentration. We assume that is the concentration for which the 
reproduction time is tr. Therefore: 

E,N = 1/tr (42) 

and Eq. ^ is satisfied. 

Similarly, EnN is the nutrition consumption rate per bacterium. We suggest that 
during reproduction time, a single bacterium consumes an amount of nutrition three 
times its mass, which is about 3 x lO^^^g. Therefore, the rate of nutrition consumption 
is: 

^ 3 X 10"^^ g nutrition 1 

EnN — (43) 

bacteria 25min 



From (^,^) we obtain the bacterial concentration: 

B = N/TEnN = 3bacteria/ij,'^ = 3 x 10%acteria / X"^ (44) 



The discrete time step of the numerical integration is measured in units of T. Sometimes 
numeric stability demands that the time step will be less than one. Then, as an example, a 
time step of 0.001 will correspond to O.OOIT ~ Isec 



C. Equations with cutoff in the reaction term 

The scaling performed in the previous section is also suitable for the cutoff equations 
(section [IV B| ). We are interested in the meaning of the cutoff (3. Since (3 is measured 



in units of B, we use (^41) to translate numerical values to bacteria concentrations. For 
example, f3 = 3x 10^^ corresponds to Ibacteria / X"^ , while (3=1 corresponds approximately 
to bacteria covering the agar surface in a continuous layer. 
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D. Mimura et al. model 



The model's equations in dimensional form (section [IV E|) are: 

I = D,V% + E,hn - . /f. (45) 
ot {an + n){ab + b) 

— = a^v^n - EM 

ot 

ds ttob 



dt {an + n){ab + b) 

where a„ and are constants which have the same units as n and b respectively. We 
introduce dimensionless variables as previously (Eq. ^3]), and impose relation ([361). We 
obtain: 

- = Db/D^V b + ^bn - ^^^j^^^^^,jB + b) 



-— = V^n — bn 



dn 
'dt 

ds aoT 



dt BN {an/N + n) {ab/B + b) 
We define new parameters: 

D = Db/D^ (47) 

We set an/N = 1 and ab/B = 1 in the last term of the bacterial equation. (The assignment 
is acceptable since Mimura et al. does not justify the exact form of that term, rather they 
state that it is only its general properties that matters.) Changing to the new parameters 
gives the dimensionless equations. 

Scaling this model requires adjusting the values of the units. Relation (|37D is replaced 
by relation (^), which is equivalent to: 

TNEb = e (50) 
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Since (^21) is still valid, we leave unchanged. The other units re-scale according to: 

T ^ eT (51) 
X ^ v^X 

B B/t 

compared to the values of the units evaluated in section ([VI B|) 
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FIGURES 

FIG. 1. Typical example of branching growth of the T morphotype for Ig/l peptone level and 
1.5% agar concentration. 

FIG. 2. Examples of typical patterns of T morphotype for intermediate agar concentration 
(clockwise from top left), (a) At very high peptone level (peptone 12g/l, agar concentration 1.75%) 
the pattern is compact, (b) At high peptone level {3g/l, agar 2%) the pattern is of dense fingers with 
pronounced radial symmetry - similar to patterns observed in Hele-Show cell, (c) At intermediate 
peptone level {Ig/l, agar 1.75%) the pattern is "bushy" fractal-like pattern, with branch width 
smaller than the distance between branches, (d) At low peptone level {O.lg/l, agar 1.75%) there 
are fine radial branches with apparent circular envelope. 

FIG. 3. (left) Fractal pattern for O.Olg/l peptone level and 1.75% agar concentration, (right) 
Dense branching pattern for 4g/l peptone and 2.5% agar. Note that the branches are much thinner 
than those in Fig. 2b, i.e. the branches are thinner for higher agar concentrations. 

FIG. 4. Patterns of colonies of T morphotype (clockwise from top left), (a) Pattern of concen- 
tric rings superimposed on a branched colony for 2.5g/l peptone level and 2.5% agar concentration. 

(b) Concentric rings in a compact growth for 15g/l peptone level and 2.25% agar concentration. 

(c) Weak chirality (global twist of the branches) for 4g/l peptone and 2.5% agar concentration. 

(d) Closer look at the branches (x50 magnification) shows density variations within each branch. 
Darker colors represent thicker layer of bacteria. 

FIG. 5. Coexistence of two morphologies near the critical point. The colony shows combination 
of the radial symmetry morphology and the fractal-like morphology. The colony is grown on agar 
concentration of 1.75% and Ig/l peptone level. 

FIG. 6. Morphology diagram of B. subtilis colonies grown by Matsushita et al. |23,33,26|. 
Taken with permission from |6C]. 
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FIG. 7. Colony of B. subtilis. Pattern of concentric rings superimposed on a branched colony. 



Taken with permission from |6C] 



FIG. 8. Electron microscope observation of T bacteria. Round or oval shapes with bright 
center are spores. Elongated shapes are living cells. The cells engolfing oval shapes are pre-spores. 

FIG. 9. Reaction term f{b) (above) and Landau-Ginzburg free energy <i>(6) (below) for the 
Fisher-Kolmogorove equation. 

FIG. 10. Typical front obtained for the ID diffusive Fisher-Kolmogorove model. Parameters: 
Db = 0.1, Dn = 1, Initial food concentration tiq = 1. 

FIG. 11. Landau-Ginzburg free energy $(6) for the solid-liquid phase- field model. Note the meta-stable 
point at 6 = and the stable one at 6 = 1. 

FIG. 12. Growth term f{b) (above) and Landau-Ginzburg free energy $(5) (below) for the 
Kessler-Levine correction. Cutoff value /3 = 1, nutrient level n = 1. 

FIG. 13. 2D growth pattern of the Kessler-Levine correction with no death term. Parameters are: 
m = 0.01, /3 = 0.25, no = 1 

FIG. 14. 2D growth pattern {b + s) of the Kessler-Levine correction with a death term. Parameters 
are: Db = 0.01, /3 = 0.25, uq = 1, /i = 0.01. 

FIG. 15. 2D growth patterns (b + s) of the model with lubrication. a. Parameters are: 
7 = -1/2, A = 0.01, r; ^ 2,Db ^ 0.5, A = 0.5, /i = l,no = 1.5, j = 0.01,;„a:r = 2, F = 1 b. Parame- 
ters are as in Fig. a, except a higher absorption rate A — 0.1. c. Parameters are as in Fig. a, except 
a different exponent for the bacterial diffusion 7 = — 1. d. Parameters are as in Fig. c, except a higher 
absorption rate A = 0.03 

FIG. 16. Profile of the fronts of the bacterial field (solid line) and of the lubricating field (dashed line) 
from the 2D model with lubrication. The bacterial field was scaled by a factor of 40. left: Parameters are 
as in Fig. [l5|a. right: Parameters are as in Fig. ^5|c. 
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FIG. 17. ID front obtained for the Kitsunezaki model. Parameters are: Dq = 0.1, k = = 0.15, tiq = 1 

FIG. 18. ID front obtained for the Kitsunezaki model with k—2. All other parameters as in Fig. 

FIG. 19. 2D growth pattern (6 + s) of the Kitsunezaki model. Parameters are: 
Do = 0.1, fc = 1,^ 0.15, no = 1. 

FIG. 20. 2D growth pattern (6) of the Kawasaki et al. model. Parameters are: Df, ^ n * b,nQ = 0.71. 

FIG. 21. Growth term f{b) (above) and Landau-Ginzburg free energy $(6) (below) for the Mimura et 
al. model. Parameters are: e = 20, /i — 2400, n = 5. Note the stable point at 6 = 0. 

FIG. 22. ID front obtained for the Mimura et al. model. Parameters are: 

Db ^ 0.1, no = 10, e = 20, fi = 2400. 

FIG. 23. Various 2D patterns {b + s) obtained for the Mimura et al. model: a. DLA-like 
{Db — 0.05, no — 10). b. dense branches {Di, — 0.09, no = 10). b. Goncentric rings (Db = 0.05, no — 12). c. 
Gompact growth {Db — 0.1, no = 14). In all cases, /i = 2400, e — 20. 

FIG. 24. Mimura et al. model: The effect of changing the initial nutrient level no. For all pictures 
Db = 0.05, = 2400, e = 20. 

FIG. 25. Mimura et al. model: The effect of changing Di,. For all pictures no = 10, ^ = 2400, e = 20. 

FIG. 26. Golonial growth velocity vs. Db for the Mimura et al. model. Parameter values as in Fig. 

FIG. 27. 2D growth pattern {b + s) of the Kitsunezaki model with food chemotaxis included, xof = 3. 
Other parameters as in Fig. ^ 

FIG. 28. 2D growth pattern {b + s) of the Kitsunezaki model with repulsive chemotactic signaling 
included. Parameters are: Xor = 1; -Dr = 1, Pr = 0.25, fl,. = 0, = 0.001. Other parameters as in Fig. |l^. 



50 



FIG. 29. 2D growth pattern {b + s) of the Mimura et al. model with food chemotaxis included. 
Xof — 0.06. Other parameters as in Fig. a. 

FIG. 30. 2D growth pattern {h + s) of the Mimura et al. model with repulsive chemotactic signaling 



included. Parameters are: xor ~ 0.1, Dr = 1, = 0.2, Q.r = 0, A,. = 0.01. Other parameters as in Fig. 23 a 
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